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Collisional  energy  transfer  modeling 
in  non-equilibrium  condensing  flows 
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Abstract.  Molecular  dynamics  and  DSMC  models  for  energy  transfer  in  argon  dimer  -  argon  monomer  collisions  are 
compared.  The  inelastic  collision  number  for  these  collisions  is  evaluated,  and  found  to  decrease  significantly  when  gas 
temperature  increases. 

Keywords:  Homogeneous  condensation,  energy  transfer,  molecular  dynamics,  DSMC 
PACS:  51.10.+y 


INTRODUCTION 

Theoretical,  experimental,  and  computational  studies  of  homogeneous  condensation  has  drawn  much  attention  over 
the  last  decades,  mainly  because  this  phenomenon  plays  an  important  role  in  many  atmospheric  and  technological 
processes,  and  understanding  of  its  physical  mechanisms  and  dependencies  is  critical  for  a  number  of  engineering 
applications.  One  of  such  applications,  pertaining  to  post  boost  vehicle  operations  at  very  high  altitudes,  is  related  to 
thruster  plume  expansion  into  surrounding  rarefied  atmosphere  [1],  Condensation  in  the  rapidly  expanding  flows  has 
been  observed  experimentally  as  early  as  mid-30s  [2],  and  has  been  extensively  studied  in  the  following  decades  (see 
for  example  [3]  and  the  references  therein).  Computational  modeling  of  expanding  condensing  flows  has  a  shorter, 
although  still  a  respectable  history.  In  the  past,  two  different  approaches  have  been  used  to  describe  homogeneous 
condensation  and,  in  particular,  cluster  nucleation  (formation  of  small  clusters  from  monomers)  in  non-equilibrium 
environment  of  rapid  expansions.  The  first  approach  is  based  on  the  classical  nucleation  theory,  and  is  in  fact  a 
macroscopic  approach  that  uses  equilibrium  thermodynamics  and  the  principle  of  detailed  balance  to  calculate  the 
cluster  size  evolution  and  properties. 

The  second  approach  is  a  kinetic,  microscopic  approach,  often  based  on  the  direct  simulation  Monte  Carlo  (DSMC 
method)  to  describe  the  flow  behavior.  The  DSMC  method  has  been  used  to  study  the  process  of  cluster  formation  and 
evolution  for  a  number  of  years.  Over  the  last  decade,  it  has  been  extensively  and  successfully  applied  to  modeling  the 
processes  of  cluster  formation  and  evolution  in  supersonic  jets  by  Levin  et  al.  (see,  for  example,  [4,  5]).  The  model 
initially  was  based  on  the  classical  nucleation  theory,  with  the  new  clusters  being  formed  at  the  critical  size.  Further 
work  of  these  authors  [4]  extended  the  kinetic  dimer  formation  approach  of  Ref.  [6],  where  it  was  assumed  that  a 
ternary  collision  always  results  in  a  dimer  formation,  to  include  molecular  dynamic  (MD)  simulations  for  obtaining 
information  on  the  probability  of  dimer  formation  in  such  ternary  collisions. 

More  recently,  the  first  principle  kinetic  theory  was  used  [7]  to  construct  a  DSMC  based  model  that  uses  a  truly 
kinetic  RRK  algorithm  [8]  to  characterize  the  cluster  evaporation  rates,  and  an  energy  dependent  collision  procedure 
similar  to  the  recombination  reaction  model  of  Ref.  [9]  for  the  collision  complex  formation.  An  empirical  parameter 
was  introduced  for  the  inelastic  collision  number  in  the  cluster-monomer  collisions,  and  the  Larsen-Borgnakke  (LB) 
principle  was  extended  to  simulate  the  energy  transfer  in  collisions  among  monomers  and  clusters.  This  principle  was 
also  used  for  the  energy  redistribution  in  the  dimer  formation  process.  The  model  [7]  was  then  adapted  and  extended  to 
a  combined  Lagrangian-Eulerian  approach  [10],  that  maintain  the  benefits  of  a  kinetic  treatment  of  the  homogeneous 
nucleation  process  while  having  a  significantly  higher  computational  efficiency  than  a  DSMC  based  method. 

The  recent  studies  have  shown  that  while  the  application  of  the  kinetic  approach  has  a  number  of  benefits,  which 
include  but  limited  to  capturing  non-equilibrium  features  of  the  condensation  process  and  correct  prediction  of  many 
condensation  phenomena  and  rates  (see,  for  example.  Refs.  [7,  11]),  they  do  suffer  from  large  uncertainties  in  such  key 
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parameters  that  they  use  as  heat  capacities,  binding  energies,  energy  redistribution  mechanism,  and  cluster  inelastic 
collision  number.  Some  physical  parameter  uncertainties  have  been  clarified  in  the  previous  paper  [11],  although  the 
energy  redistribution  mechanisms  in  cluster  formation  and  monomer-cluster  collisions  are  still  largely  unknown. 

The  main  objective  of  this  work  is  the  study  of  the  energy  transfer  between  translational  and  internal  modes  in 
collisions  between  small  argon  clusters  and  monomers,  evaluation  of  the  applicability  of  the  Larsen-Borgnakke  model 
to  simulate  this  transfer,  and  analysis  of  the  impact  of  this  transfer  on  expanding  plume  flows.  The  molecular  dynamics 
method  is  used  to  compute  energy  redistribution  after  argon  cluster  -  monomer  collisions,  and  the  Lagrangian-Eulerian 
approach  [10]  is  employed  for  the  plume  flow  analysis. 


IMPACT  OF  THE  INTERNAL  ENERGY  TRANSFER  RATE 

It  has  been  previously  noted  [7]  that  the  transfer  of  energy  between  translational  and  internal  modes  in  collisions 
between  dimers  and  monomers  has  a  strong  impact  on  the  dimer  equilibrium  constant.  The  reason  for  this  effect  is 
that  the  dimer-monomer  collisions  are  the  principal  factor  that  changes  dimer  internal  energy.  Dimer  internal  energy 
distribution,  in  turn,  directly  affects  the  dimer  dissociation  rate  and,  therefore,  the  equilibrium  constant.  The  dimer 
dissociation  is  an  endothermic  process,  and  may  be  considered  as  a  result  of  unimolecular  dissociation  of  At'2  clusters 
that  have  internal  energy  larger  than  the  binding  energy  of  the  clusters.  It  is  obvious  that  the  lower  the  dimer  binding 
energy,  the  larger  the  impact  of  the  dimer  internal  energy  relaxation  rate.  It  was  found  in  Ref.  [7],  where  the  cluster 
internal  energy  relaxation  number  Z  was  introduced,  similar  to  the  rotational  and  vibrational  relaxation  numbers  of 
gas  molecules,  that  Z  has  a  much  stronger  influence  on  the  argon  cluster  nucleation  process  than  on  water  cluster 
nucleation.  This  was  related  to  the  Ar2  binding  energy  Ej,  of  1.98  x  10  21  J,  as  compared  to  water  dimer  Ej,  of 
1.98  x  10  21  J.  Even  though  water  dimers  have  a  significantly  larger  heat  capacity  than  argon  dimers,  it  still  takes  on 
average  much  more  collisions  to  increase  the  internal  energy  Ejnt  of  a  stable  water  dimer  that  is  initially  formed  with 

E in!  E  Ely. 


FIGURE  1.  Impact  of  cluster  internal  energy  relaxation  rate  on  equilibrium  constant  and  comparison  with  theoretical  predictions 
[12]. 

For  stable  argon  dimers,  it  typically  takes  only  few  collisions  to  exceed  Ej,.  Although  dimers  with  after-collision 
internal  energy  larger  than  the  binding  energy  do  not  dissociate  immediately  after  collision,  their  lifetime  (the  RRK 
model  [13]  may  be  used  to  estimate  the  cluster  lifetime)  is  typically  shorter  than  the  mean  collision  time.  This  makes 
the  dissociation  of  dimers  with  Elnt  >  Ej,  likely.  The  impact  of  Z  on  the  argon  dimer  equilibrium  constant  Keq  is 
illustrated  in  Fig.  1,  where  the  model  [7]  is  used  to  calculate  Keq  in  a  thermal  bath  of  pure  argon  with  a  density 
of  1.326  kg/m3.  Note  that  the  dependence  on  gas  density  in  this  case  is  negligible,  since  the  mean  collision  time  is 
much  larger  than  the  average  cluster  lifetime.  The  value  of  Z  was  varied  from  1  to  12.5  for  two  gas  temperatures.  As 
expected,  the  increase  of  Z  results  in  significant  increase  of  the  equilibrium  constant.  The  dependence  on  Z  is  almost 
linear  at  small  Z,  but  becomes  weaker  as  Z  increases.  This  is  because  the  monomer-dimer  sticking  becomes  dominant 
for  Z  1  (the  sticking  coefficient  was  0.06  for  this  case).  As  seen  from  the  figure,  the  selection  of  an  optimum  value 
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TABLE  1.  Inelastic  collision  number  for  different  gas  tem¬ 
peratures. 


T.  K 

0.0 

100.0 

200.0 

300.0 

400.0 

500.0 

Z'1 

0.25 

0.13 

0.08 

0.06 

0.046 

0.04 

of  Z  that  provides  best  agreement  with  the  theoretical  results  of  Ref.  [12]  depends  on  thermal  bath  temperature.  It  is 
about  7  for  100  K  and  12  for  200  K.  The  temperature  dependence  used  in  Ref.  [11]  is  summarized  in  Table  1. 


NUMERICAL  APPROACH 

A  general  chemical  dynamics  program  VENUS  [14]  has  been  used  in  molecular  dynamics  simulations.  The  Morse 
pairwise  potential  was  applied  to  described  the  inter-  and  intra-cluster  interactions  of  argon  atoms.  Note  that  although 
this  is  a  less  accurate  interaction  potential  than  a  more  sophisticated  HFDID1  potential  [15],  it  has  very  similar  bonding 
energy,  and  is  expected  to  provide  qualitatively  similar  results.  The  Morse  potential  is  represented  in  its  standard  form 

V(r)  =  D(exp(-2a(r  —  re))  — 2exp  (—a(r  —  re)) 

with  the  following  parameters  taken  from  [16]:  a=0. 908597,  D= 99  cm-1,  rf=3.757  A.  Trajectories  started  at  a 
separation  of  13Abetween  dimer  center  of  mass  and  impinging  atoms.  The  convergence  of  cross  section  for  various 
impact  parameters  and  number  of  trajectories  was  studied,  and  for  typical  runs  maximum  impact  parameter  of  lOAwas 
used,  and  20000  trajectories  were  run.  Trajectories  were  terminated  when  the  distance  between  the  dimer  center  of 
mass  and  impinging  atom  was  13  A. 

The  gas  flow  was  simulated  with  the  DSMC  code  SMILE  [17]  where  the  first-principle  model  [7,  11]  was  incorpo¬ 
rated.  Some  changes  were  introduced  in  the  parameters  of  the  approach  that  are  believed  to  more  accurately  reproduce 
the  physics  of  the  dimer  formation  phenomenon.  Namely,  Lennard-Jones  parameters  used  in  the  expression  [18]  to 
estimate  the  collision  complex  lifetime  are  e  =  1 .98  x  1 0  2 1  J  and  a  =  3.405  A.  The  monomer  vibration  frequency  in 
argon  clusters  was  found  to  be  2.6  x  1011  s-1.  Finally,  the  dimer  stabilization  probability  of  0.2  was  used. 

AFTER-COLLISION  ENERGY  REDISTRIBUTION 

Previous  model  [7]  assumed  the  Larsen-Borgnakke  type  energy  redistribution  after  inelastic  argon  monomer-dimer 
collisions,  with  the  inelastic  collision  number  Z  fitted  to  match  the  dimer  formation  equilibrium  constant  as  a  function 
of  temperature  (see  Table  1).  While  this  approach  appears  reasonable  for  preserving  correct  dimer  mole  fractions  at 
equilibrium,  it  does  not  guarantee  proper  energy  distributions,  trimer  and  larger  cluster  formation  rates,  or  even  dimer 
formation  rates  in  non-equilibrium  plume  flows.  In  order  to  analyze  the  reliability  and  applicability  of  this  approach, 
more  accurate  analysis  is  necessary,  such  as  that  of  trajectory  calculations  with  a  realistic  interaction  potential. 

The  profiles  of  after-collision  internal  (the  sum  of  rotational  and  vibrational)  energy  of  argon  dimers  is  presented 
in  Fig.  2  for  two  different  approaches,  the  molecular  dynamics  trajectory  calculations  and  the  Larsen-Borgnakke 
collision  algorithm.  The  calculations  were  conducted  for  the  initial  vibrational  level  v  =  2  and  rotational  level  /  =  20 
and  the  Maxwellian  distribution  of  relative  collision  velocities  at  a  temperature  of  T *  =  2.  Here  T*  is  the  temperature 
normalized  by  the  reduced  argon  dimer  bond  energy  of  E^/k  =  143.2  K.  Only  dimer-monomer  collisions  that  resulted 
in  net  dimer  internal  energy  change  were  considered  for  this  plot.  It  is  seen  that  there  is  a  qualitative  difference  in 
the  after-collision  distributions,  as  the  MD  collision  tend  to  cause  only  a  small  change  in  the  dimer  internal  energy, 
whereas  the  LB  approach  results  in  a  Boltzmann  distribution  at  a  local  temperature  of  about  185  K. 

Note  that  stable  dimers  are  those  with  the  internal  energy  smaller  that  Eb,  that  is  the  particles  that  fall  to  the  left 
of  Eb/k  =  143.2  K  in  the  after-collision  energy  distribution  profile.  The  dimers  with  internal  energy  in  excess  of  Ef, 
are  unstable  and  may  dissociate  before  the  next  binary  collision.  Note  that  the  time  to  dissociation  may  be  estimated, 
for  example,  using  the  RRK  model  or  more  sophisticated,  and  is  generally  on  the  order  of  the  inverse  dimer  vibration 
frequency  of  2.6  x  1011  s-1.  For  a  1  atm  and  300  K  argon  gas,  where  the  dimer  mean  collision  time  is  on  the  order 
of  10“10  S,  most  of  unstable  dimers  will  dissociate  before  the  next  collision.  In  order  to  create  the  same  number  of 
unstable  dimers  after  a  single  collision,  a  parameter  Z  needs  to  be  used  in  the  LB  procedure,  so  only  a  fraction  1/Z 
leads  to  the  internal  energy  transfer.  The  use  of  Z  will  result  in  a  proportional  decrease  of  the  dashed  line  in  Fig.  2,  but 
will  not  change  the  equilibrium  shape  of  that  curve,  characterized  by  relatively  large  population  of  high-energy  tail. 
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FIGURE  2.  After-collision  internal  energy  redistribution  in  inelastic  LB  and  MD  collisions.  Initial  internal  energy  levels  are  v=2, 
J=20  (E /k  =  83.4  K. 


DIMER  BREAKUP  CROSS  SECTIONS  IN  DIMER-MONOMER  COLLISIONS 

The  argon  dimers  are  characterized  by  fairly  low  binding  energy,  which  results  in  high  probability  of  formation  of 
unstable  dimers  after  dimer-monomer  collisions,  followed  by  dimer  breakup.  Therefore,  it  is  important  to  analyze 
the  stable-to-unstable  (or  dimer  breakup)  collision  cross  sections  at  different  gas  temperatures.  From  these  collision 
cross  sections,  obtained  with  the  MD  approach,  the  inelastic  collision  number  for  dimer-monomer  collisions  may  be 
evaluated  for  a  given  intermolecular  potential  that  specifies  the  total  collision  cross  section  (such  as  the  VHS  model  in 
DSMC).  Consider  first  the  impact  of  the  vibrational  level  v  and  the  rotational  level  /  of  a  colliding  argon  dimer  on  the 
cross  section  of  a  process  that  transforms  a  stable  dimer  with  the  internal  energy  Elnt  <  to  an  unstable  dimer  with 
Ejn,  <  E],.  The  MD  results  obtained  for  different  v  and  J  are  shown  in  Fig.  3.  The  obtained  cross  sections  are  close 
to  the  total  collision  cross  section  of  80Acalculated  using  the  VHS  model.  As  expected,  the  cross  section  generally 
increases  with  rotational  and  vibrational  level,  which  reflects  the  increase  of  the  dimer  breakup  probability  with  the 
available  internal  energy  of  the  colliding  dimer.  Note  that  the  stable-to-unstable  dimer  cross  section  is  more  sensitive  to 
the  change  in  vibrational  level.  Fig.  3  (left),  where  the  energy  quantum  is  larger;  the  cross  section  increase  is  relatively 
insignificant  for  small  rotational  numbers.  Fig.  3  (right),  since  the  J 2  dependence  of  the  internal  energy. 
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FIGURE  3.  Stable-to-unstable  dimer  cross  section  as  a  function  of  changing  vibrational  (left)  or  rotational  (right)  level.  T*  =  0.5. 

The  increase  of  the  cross  section  with  the  rotational  quantum  number  is  more  pronounced  at  lower  temperatures,  as 
illustrated  in  Fig.  4  (left).  This  is  related  to  the  fact  that  at  T *  =  2,  the  average  relative  translational  energy  is  noticeably 
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larger  than  Eb,  and  thus  the  pre-collision  rotational  energy  is  relatively  less  important  than  for  T*  =  0.5.  The  average 
internal  energy  of  stable  dimers  changes  between  Ejnt/k  =  78  K  and  Ejnt/k  =  92  K,  and  the  dimer  breakup  cross 
section  for  a  typical  internal  energy  of  84  K  (v  =  2  and  J  =  20)  is  shown  in  Fig.  4  (right).  It  is  interesting  to  note 
that  the  breakup  cross  section  practically  does  not  change  after  T*  «  1.  This  indicates  that  dimer  breakup  probability, 
which  is  the  ratio  of  the  breakup  cross  section  to  the  total  collision  cross  section  G,ot ,  will  only  be  a  function  of  atot 
at  T  «  140 K.  For  a  lower  internal  energy  (v  =  2  and  J  =  4),  the  temperature  dependence  is  more  noticeable,  but  still 
relatively  weak. 
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FIGURE  4.  Stable-to-unstable  dimer  cross  section  for  v  =  2  and  two  different  translational  temperatures  (left)  and  for  two 
different  internal  states  and  varying  translational  temperature  (right). 


DIMER  BREAKUP  AND  INELASTIC  PROBABILITIES 

The  most  important  outcome  of  collisions  between  stable  droplets  and  monomers  is  the  formation  of  unstable  dimers 
with  internal  energy  larger  than  the  binding  energy,  since  such  a  formation  is  the  main  mechanism  of  droplet  breakup. 
Its  impact  on  dimer  population  will  in  turn  influence  the  population  of  larger  clusters,  that  are  created  through  sticking 
and  coalescence  of  dimers.  For  the  condensing  gas  modeling  to  be  reliable,  it  is  important  to  implement  correct  rates 
of  droplet  breakup.  If  one  assumes  that  for  maintaining  correct  rates  of  droplet  breakup  it  is  sufficient  to  have  realistic 
number  of  stable-to-unstable  dimer  collisions  for  any  given  internal  energy,  then  the  key  issue  is  to  maintain  correct 
droplet  breakup  cross  sections. 

For  a  kinetic  condensation  model  that  is  based  on  the  VHS  total  collision  cross  sections  and  the  LB  energy 
redistribution  algorithm,  that  means  that  the  value  of  the  inelastic  collision  number  Z  should  be  chosen  from  the 
relation  Z~l  =  (j)  =  Or(,f/OiiS.  where  0  is  the  probability  of  the  collisional  energy  transfer  that  results  in  the  formation 
of  an  unstable  dimer,  Gref  is  the  reference  dimer  breakup  cross  section  that  may  be  obtained,  for  example,  from 
molecular  dynamics  simulations,  and  Gib  is  the  dimer  breakup  cross  section  in  a  single  VHS  collision  that  has  an 
LB  energy  transfer  between  all  translational  and  internal  energy  modes.  The  results  of  calculation  of  (!)  for  different 
rotational  and  vibrational  levels  and  translational  temperatures  are  summarized  in  Table  2. 

Two  conclusions  may  be  drawn  from  that  table.  First,  the  MD  calculations  indicate  that  for  a  given  internal  energy, 
the  dimer  with  a  higher  vibrational  energy  will  break  up  with  higher  probability  (vibrational  favoring).  Second,  and 
most  importantly,  for  given  vibrational  and  rotational  levels,  the  inelastic  collision  number  Z  defined  as  Z  =  Gm / Gmd 
significantly  increases  when  temperature  increases.  For  example,  for  the  reduced  internal  energy  of  83.89  K,  (j)  =  1/Z 
decreases  by  almost  a  factor  of  4,  from  1.18  to  0.3,  when  the  gas  temperature  increases  from  T*  =  0.5  to  2.  Similar 
temperature  dependence  is  observed  when  Z  is  calculated  in  a  conventional  manner  based  on  energy  relaxation  to 
equilibrium  value.  Note  that  for  larger  dimer  internal  energies,  the  values  of  (j)  become  larger  than  unity,  which  is 
related  to  incorrect  after-collision  energy  redistribution  in  the  LB  approach,  as  was  discussed  in  the  previous  section. 

While  the  values  of  (j)  calculated  as  (f>  =  Gmd/Glb  may  be  directly  used  in  modeling  condensing  flows  with  a  kinetic 
approach  such  as  that  of  Ref.  [7],  the  use  of  the  LB  technique  for  energy  transfer  in  dimer-monomer  collisions  results 
in  drastically  different  from  MD  energy  redistribution,  which  may  lead  to  inaccurate  results.  Therefore,  it  is  better 
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TABLE  2.  Droplet  breakup  cross  sections  for  different  internal 
energies 


V 

J 

¥.K 

&MD,  A2 

Gvhs,  A1 

d LB,  A2 

0 

0 

0 

11.71 

r*  = 

9.08 

=  0.5 

112.64 

15.22 

0.60 

0 

10 

20.87 

10.71 

112.64 

16.70 

0.64 

0 

20 

46.11 

15.68 

112.64 

21.45 

0.73 

0 

30 

85.60 

32.83 

112.64 

30.92 

1.06 

2 

0 

53.53 

20.55 

112.64 

23.03 

0.89 

2 

10 

61.64 

23.31 

112.64 

24.87 

0.94 

2 

20 

83.89 

36.06 

112.64 

30.45 

1.18 

2 

30 

115.82 

70.07 

112.64 

39.79 

1.76 

4 

0 

87.32 

43.29 

112.64 

31.39 

1.38 

4 

10 

94.32 

48.69 

112.64 

33.37 

1.46 

4 

20 

113.25 

71.88 

112.64 

38.99 

1.84 

2 

0 

53.53 

T * 
37.64 

=  2 

226.26 

160.88 

0.23 

2 

10 

61.64 

40.36 

226.26 

162.80 

0.25 

2 

20 

83.89 

50.13 

226.26 

167.99 

0.30 

2 

30 

115.82 

81.78 

226.26 

175.06 

0.47 

to  evaluate  (f>  from  comparison  with  equilibrium  dimer  formation  rates,  and  select  the  values  that  provide  a  good 
agreement  with  available  theoretical  results.  It  was  found  from  comparison  with  theoretical  analysis  [12]  that  0  =  0.2 
for  equilibrium  gas  at  T*  =  0.7,  0.11  for  T*  =  1.4,  and  0.08  for  T*  =  2.1.  This  presents  a  temperature  dependence 
similar  to  that  obtained  with  0  =  Omd/Glb*  although  smaller  magnitude  of  the  inelastic  collision  probability. 


CONCLUSIONS 

The  molecular  dynamics  calculations  of  energy  transfer  in  argon  dimer-monomer  collisions,  with  emphasis  on  col¬ 
lisions  that  result  in  formation  of  unstable  dimers  with  internal  energy  larger  than  the  binding  energy.  Comparison 
with  the  energy  transfer  algorithm  based  on  VHS  and  Larsen-Borgnakke  models  for  total  and  inelastic  cross  sections 
indicate  that  the  inelastic  collision  number  decreases  strongly  with  temperature.  The  large  difference  between  the 
after-collision  internal  energy  distributions  in  MD  and  LB  models,  as  well  as  the  presence  of  vibrational  favoring  in 
MD,  indicate  to  potential  problems  when  using  the  LB  model  for  energy  transfer  in  droplet-monomer  collisions. 
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